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The problem of determining the coefficients for a linear feedback 
guidance system with polynomial gains is investigated. The solution to 
the linearized equations of motion is employed to determine the first 
partial derivatives of the optimum thrust angle deviations with respect 
to position and velocity coordinates. The partial derivatives are 
approximated by quadratic time functions, which in turn are the first 
approximation to a time-variable set of gains. This system is then 
analyzed, and two of the weighting functions are linearized with respect 
to perturbations in the polynomial coefficients. Changes in these coef- 
ficients are determined to minimize the weighted sum of squares of these 
expressions over several time points. The resulting changes were made, 
and a satisfactory reduction in magnitude of these weighting functions 
resulted. This, coupled with a linear feedback of thrust acceleration 
and a small time function derived to cancel the effect of initial con- 
ditions, produced the final guidance function which performed exception- 
ally well on a 100 n.m. orbital mission of an early SA-6 second stage 
vehicle. Of particular interest is a demonstration of the effectiveness 
of the numerical methods used which are generally applicable to a wide 
variety of guidance and control feedback problems. 
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LINEAR FEEDBACK GUIDANCE 
SUMMARY 


The problem of determining the coefficients for a linear feedback 
guidance system with polynomial gains is investigated. The solution to 
the linearized equations of motion is employed to determine the first 
partial derivatives of the optimum thrust angle deviations with respect 
to position and velocity coordinates. The partial derivatives are 
approximated by quadratic time functions, which in turn are the first 
approximation to a time-variable set of gains. This system is then 
analyzed'', and two of the weighting functions are linearized with respect 
to perturbations in the polynomial coefficients. Changes in these coef- 
ficients are determined to minimize the weighted sum of squares of these 
expressions over several time points. The resulting changes were made, 
and a satisfactory reduction in magnitude of these weighting functions 
resulted. This, coupled with a linear feedback of thrust acceleration 
and a small time function derived to cancel the effect of initial con- 
ditions, produced the final guidance function which performed exception- 
ally well on a 100 n.m. orbital mission of an early SA-6 second stage 
vehicle. Of particular interest is a demonstration of the effectiveness 
of the numerical methods used which are generally applicable to a wide 
variety of guidance and control feedback problems. 


INTRODUCTION 


The problem of determining a closed loop guidance function with 
linear, time-variable feedback is investigated by studying the behavior 
of the equations of motion in the neighborhood of a standard calculus of 
variations solution. The linearized Euler-Lagrange equations are used 
to determine a good approximation to the optimum thrust angle deviation, 
AX. The solution gives AX as a linear combination of known functions of 
time. The particular linear combination required is that one which ful- 
fills the desired end conditions. The linearly predicted deviations of 
the end conditions, as a function of deviations in initial conditions, 
are used to determine the first partial derivatives of the optimum thrust 
angle with respect to position and velocity components. Quadratic time 
approximations to these partial derivatives provide the time -variable 
gains for a linear feedback system which is analyzed by the same methods 
employed in the open loop analysis described ip Reference 1. To reduce 
the effect that thrust acceleration variations would have on this system. 



the weighting functions involved are linearized with respect to changes 
in the coefficients of the gain polynomials. Changes in these coeffi- 
cients are determined to minimize the weighted sum of squares of the 
resulting expressions over a number of time points. Although this is 
only one step of an iterative procedure, the reduction in magnitude of 
the weighting functions from this one iteration is sufficient to demon- 
strate the effectiveness of the method. What little effect of thrust 
acceleration deviations is left is further reduced by the addition of 
thrust acceleration feedback with gains which are also quadratic time 
functions. A polynomial of time, whose coefficients are determined as 
linear functions of initial conditions, effectively cancels the remaining 
error due to initial conditions. The resulting guidance function has been 
determined for an early SA-6 second stage vehicle designed for a 100 n.m. 
circular orbit. The tabulated results for a number of examples reflect 
the effectiveness with which linear analysis can be applied to problems 
of this type. 


I. LINEARIZED EULER-LAGRANGE EQUATIONS 


An explicit solution to the linearized equations of motion was 
determined in Reference 1. The deviations in end conditions, Ar and A0, 
and the additional burning time, At, were determined as functions of 
deviations in initial conditions, AX,-, , variations in thrust angle, AX, 
and variations in thrust acceleration, Af/m. These terminal deviations 
were found to be adequately represented by the following expressions: 


LI 

At=U 0 (t 0 ) AX(t 0 ) + J h(AX, Af/m, t) dt. 


( 1 . 1 ) 


ii 

Ar = U!(t 0 ) AX(t c ) + J f (AX, Af/m, t) dt, 


( 1 . 2 ) 


-n 


A9 = U 2 (t 0 ) AX(t 0 ) + J g(AX, Af/m, t) dt. 


(1.3) 


2 


where 


h(AX, M/m, t) = ho M/m + 


f (AX, M/m, t) = f G M/m + 


g(AX, M/m, t) = g Q M/m + 


+ 


+ 


+ 


Af/m\ 

~7frJ 


Af/m\ 

f/mj 


M/m\ 

f/tj 


(h x AX + h 2 AX 2 + 


(f x AX + f 2 AX 2 + ...), 


( gl AX + g 2 AX 2 + ...), 


(1.4) 


(1.5) 


( 1 . 6 ) 


Using the Euler-Lagrange equation to minimize At under the constraint 
that At and A9 have specified values leads to the following necessary 
condition: 


‘ M/m 


f/mj 



(h x + X x f x + Xegj.) + 2(h 2 + X x f 2 + >eg 2 ) AX + ...J = 0. 

(1.7) 

For the standard calculus of variations solution, this means that 


h i + 7M.fi + >egi = 0. 


( 1 . 8 ) 


To simplify further analysis. 


the following definitions are employed 



Substitution of these expressions into equation (1.7), neglecting terms 
higher than first order and taking cognizance of the relationship 
defined in equation (1.8), results in the following: 


fiAAi + giAA 2 + 2(h 2 + Aif 2 + A 2 g 2 ) AX. = 0. 


( 1 . 10 ) 


For convenience, the following quantities are defined: 


Pi = 


aLl 

loo (h 2 q* A^_f 2 q* A pg p) 


Pa = 


i&i 

(h 2 + A x f 2 q- A 2 g 2 ) 


lOOAAi 


P ~ [Pi Pal 

a = r a 

L p J 


(l.u) 


With these definitions, it can be readily verified that the follow- 
ing expression for AX satisfies equation (1.10): 


AX = PA. 


( 1 . 12 ) 


Any linear combination of the functions p, and p p will satisfy . 
to first order, the Euler-Lagrange equation . The functions p x and p 2 
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are illustrated in Figure 1.1. These particular functions can be well 
represented by the following quadratic functions of time; 


Pi = 5.8096 + .1734 t - .2948-r 2 
p 2 = -6.5992 - 3. 2636 t + .IOOOt 2 


(1.13) 


where 


T = 


t - t 

O 


100 


It was shown in Reference 1 that for this same mission the value 
of X on the standard trajectory could be adequately represented by a 
quadratic in time. Thus, for this example, X can be well represented 
as a quadratic time function with two arbitrary parameters, a and fj. 

If a and p can be determined in any way so that Ar = A0 = 0, the result- 
ing guidance function will satisfy the Euler-Lagrange equation to first 
order. 

In actual practice, a single pair of constants, a and (3, to be 
used throughout the trajectory cannot be determined because of the 
future uncertainties in the thrust vector. However, if their best 
approximation is available, a good quadratic approximation to the 
optimum X is available until such time as later information provides 
better values for a and p. As the time intervals between corrections 
in a and (3 become smaller, a and (3 become continuous; AX then becomes a 
time variable function of the current state of the system. 
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II. FIRST PARTIAL DERIVATIVES 


In order to attain a first approximation to a linear feedback sys- 
tem, equations (1.2) and (1.3) will be linearized with respect to AX, 
and all other terms will be considered zero except the initial condi- 
tions, £Xq. Then equation (1.12) can be used to determine the linear 
effect of AXq on the optimum value of AX. To simplify subsequent deriva- 
tions, the following definitions will be used: 


AR = 


Xar 

A8 


( 2 . 1 ) 


F. = 
1 



i - 0,1,2, 


( 2 . 2 ) 


M(t 0 ) 



u(t D ) 


u i( t o> 

U s( t 0 ) 


(2.3) 


(2.4) 


Considering only the effects of AXq and the linear effect of AX, 
equations (1.2) and (1.3) can be written in the following matrix form: 


AR = U(t Q ) AX(t 0 ) 



AX dt. 


(2.5) 
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Substituting the expression for AX obtained in equation (1.12) and using 
the definition of equation (2.3) yields the following expression: 

= U(t Q ) AX(t 0 ) + M(t 0 ) A . (2.6) 

Setting AR = 0 and solving for A gives 

A= -M _1 (t 0 ) U<t 0 ) AX (t Q ). (2.7) 

The expression for AX. at t = t Q then becomes 

AX(t 0 ) = -P(to) M" 1 (t 0 ) U(t 0 ) AX(t 0 ). 

The uniqueness of an existing series expansion ensures that 


-P(t 0 ) M- 1 (t () ) U(t D ) = 


dx Sx dx dx 

Sx Sy c)x ciy 


J t=t f 


and since the derivation is independent of the particular value of t D , 
the following expression is obtained for the first partial derivatives; 


c>X Sx Sx Sx 

Sc cS 1 % 


= -P(t) W i (t) U(t). 


( 2 . 8 ) 


A polynomial fit of these partial derivatives will be used as a first 
approximation to the desired linear feedback system. Subsequent analysis 
will determine the acceptability of this system, and methods of improving 
its performance will be outlined and applied. Figures 2.1 and 2.2 show 
the values obtained for the first partial derivatives from equation (2.8) 
where they are compared with their quadratic approximations defined by 
the following functions: 


not filmed. 
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a(t) 
b (t) 
c(t) 
d(t) 


= a o + ajT + a 2 x 2 ~ ^ 

= b 0 + biT + b 2 x 2 ~ 

= e 0 + CtlT + c 2 t 2 

= d Q + djT + d 2 x 2 ~|r J 


(2.9) 


a Q = .104259 

II 

H 

-.160995 

a 2 - . 07 6398 

b Q = .37021 

bi = 

-.44893 

b 2 = .21551 

c Q = .01551 

c i = 

-.003883 

c 2 = .003477 

d Q = .062320 

di = 

-.010020 

d 2 = .007364 


( 2 . 10 ) 


t - t 


T = 


100 


The deviation in thrust angle defined by this system can be expressed 
as follows: 


AX = BAX + 6X , 


( 2 . 11 ) 


where 


B = [a(t) b(t) c(t) d(t)]. 


( 2 . 12 ) 


The quantity SX represents any disagreement between the commanded thrust 
angle and that which actually exists. 

The effect on end conditions of deviations in the thrust accelera- 
tion vector is reflected by the weighting functions f i and g i shown in 
equations (1.5) and (1.6). A review of the derivation of these equa- 
tions in Reference 1 reveals the following relationships: 
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F 2 = - (n/180) 2 ^2 F 0 




and / 

F 3 = - (it/180) 2 | Fi 

J 

where 

> i “■ 0 j ^ ) • • . * 



(2.13) 


Furthermore, all even ordered vectors are constant multiples of 
F^, and odd ordered vectors are constant multiples of F x . It follows 
that all the weighting functions can be determined immediately from 
either F Q or F x . These two vectors represent the weighting functions 
associated with the linear terms of Af/m and AX, respectively, and can 
be determined from the solution to the following system: 


AX = AA X + H 0 Af/m + H x AX. 


(2.14) 


The solution is 


tn 


AR = U(t 0 )AX(t 0 ) + / F x AX dt 4- 


F Q Af/m dt. 


(2.15) 


Substitution of the expression for AX shown in equation (2.11) into 
equation (2.14) results in the following system: 


AX = (A + H X B) A X + H q Af/m + H x 5X 
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or 


AX = A AX + H 0 AE/m + Hi 8X, 


(2.16) 


where 


A = A + H X B. 


This expression is of identical form to that appearing in equation 
(2.14). The solution can be expressed as follows: 


AR = U(t Q ) AX(t 0 ) 



F Af/m dt 
o 


tn 


F x SX dt. 


(2.17) 


A comparison of U(t Q ) and U(t Q ) for the mission under consideration 
yields the following: 


u(t 0 ) 


.349481 

-.008175 


1.25909 .170980 .491843 

-.007129 -.005744 -.008092 


(2.18) 


and 


U(t Q ) 


.045516 .033588 

-.007110 -.000262 


.017544 

-.003038 


.00957 

-.000534 


(2.19) 


The coefficients above were obtained for variables expressed in 
the following units: 

Ac, Ay and Ar are expressed in kilometers. 

Ax and Ay are expressed in meters per second. 

A0 is expressed in degrees. 
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# 


Although the elements of U(t Q ) are not small enough to be ignored, 
they provide a good linear prediction of the effect of AX(t 0 ) on the end 
conditions* With this knowledge available, it should not be difficult 
to determine a time function of second stage ignition to cancel these 
effects. With this_in mind, primary consideration will be given to the 
elements of F Q and F^. The elements of these two vectors are illustrated 
in Figures 2.3 and 2.4 where they are compared with the corresponding 
elements of F 0 and F x . That the magnitude of the weighting functions 
defined by equations (2.9) - (2.12) is apparent. A better indication 
of the effect that 5X and Af /m might have on Ar and AO is reflected in 
the following values which have been determined. 



.498 km/deg 

.12 


( 2 . 20 ) 



4.135 km sec 2 /m 
.550 deg sec 2 /m 


( 2 . 21 ) 


Although the effect of 6X represented by equation (2.20) might be 
acceptable, it is apparent by equation (2.21) that the effect of varia- 
tions in f/m is not. Even variations as small as .25 m/sec 2 could 
cause as much as 1 km in Ar. The magnitude of elements of F are func- 
tions of the coefficients defined by_equations (2.10). By altering 
these coefficients, the elements of F 0 can be altered. A procedure 
to accomplish this will be outlined in the following section. 


i 
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III. PERTURBING THE POLYNOMIAL COEFFICIENTS 


The weighting functions represented by the vectors Fo 
determined from the following system of equations: 

and F x were 

Fi(t) = T U(t n ,t) R ± (t) (i = 0, 1, ...) 

, (3.1) 

• 

UCt.ti) = A(t) U(t, tl ), 

(3.2) 

where 


= i 

(3.3) 

and 


A(t) = A(t) + HxCt) B(t) . 

(3.4) 

The effect of changing the polynomial coefficients in B(t) can be 
investigated by defining another weighting function vector F 0 (t) resulting 
from another choice of coefficients defining a different function B(t) 
where 

F 0 - F 0 + AFq 

(3.5) 

B = B + AB. 

(3.6) 

A 

Then Fo is determined from the following relationships: 


F 0 (t) = T U(t n ,t) H 0 (t) 

(3.7) 

U(t,t!) = A(t) U(t,t!> 

(3.8) 

where 


U(ti,tx) = I 

(3.9) 

and 


A(t) = A(t) + H x (t) B(t). 

(3.10) 

, . viv NOT 

pa,GE ^ 

P * 
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Defining the new variables in terms of the previous ones yields the 
following definitions: 


A(t) = A(t) + AA (3.11) 

U(t,t x ) = U(t,t x ) + AUCt.tx). (3.12) 

From equations (3.10), (3.6), and (3.4), it can be seen that 

AA = HjAB. (3.13) 

Substitution of equations (3.11) and (3.12) into equation (3.8) yields 
U(t,t x ) +AU(t,t x ) = A(t) U(t,t x ) + A(t) AU(t,ti) + AA(t) U(t,t x ) 

+ AA(t) AU(t,t x ). 

The relationship defined by equation (3.2) reduces this to 

AU(t , t x ) = A(t) AU(t,t x ) + AA(t) U(t,t x ) + AA(t) AU(t,t x ). (3.14) 

Equations (3.3), (3.9), and (3.12) imply the following constraint on the 
initial conditions: 

AU(t x , t x ) = 0. (3.15) 

From equations (3.1), (3.5), and (3.12), it can be determined that 

AFo(t x ) = T AU(t n ,t x ) H 0 (t!), (3.16) 

where AU(t n ,t x ) satisfies equation (3.14). 
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An approximate solution, F* can be determined by linearizing equation 
(3.14). Then 


*£(ti) = FoCtiJ+A^Cti) 


where 


AB£( tl ) = T AU(t n>tl ) ^(ti). 


Aua,^) = A(t) AUCt.ti) +AA(t) U(t,t x ) (3.19) 


^(ti.tj.) = 0. 


(3.20) 


The solution to equation (3.19) can immediately be determined: 


^(tn’tl) = U(t n ,t!) AUCti.ti) + J U(t n ,t) AA(t) U(t,t x ) 

ti 


From equation (3.13) and (3.20), this can be reduced to 


AlKtn.tl) = / U(t n ,t) H x (t) AB(t) U(t ,t x ) dt. 


Inserting this expression into equation (3.18) gives the following equati 


A F 0 (tx) - / T U(t n ,t) H x (t)AB(t) U(t,t x ) ^(tj.) dt. 
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Equation (3.1) reduces this to 


AF*(t!> = / ^(t) AB(t) H 0 (t x ) dt. 


(3.21) 


AB(t) can be expressed in terms of changes in the polynomial 
coefficients as follows: 


AB(t) = [Aa(t) Ab(t) Ac (t) Ad(t)] 


Aa(t) = Aao + Aa^T + AagT 2 
Ab(t) = Ab 0 + Ab x t + Al-oT 2 
Ac(t) = Ac 0 + AciT + Acg T 2 
Ad(t) = Ad 0 + AdiT + AdaT 2 


(3.22) 


Equations (3.17), (3.21), and (3.22) provide an approximation to the 
vector F 0 expressed as a linear time-variable function of perturbation 
in the polynomial gain coefficients. This expression has been evaluated 
over a number of time points for the application being investigated. The 
least squares criterion was applied to the function. 


[f* 2 (ti) + 25g* 2 (t i )]. 


22 



In order to take advantage of an available computer program which could 
handle only a restricted number of coefficients, Aap, £b 0 , ^o* anc * ^o 
were set equal to zero and the remaining eight coefficients were deter- 
mined to minimize expression (3.23) with the following results. 




Aai = 

-.007333 

Aa 2 = 

.046081 

Ab x = 

.764653 

Ab 2 = 

-.093713 

II 

H 

* 

.045045 

A c 2 = 

-.007476 

Ad x = 

.077065 

Ad 2 = 

-.007011 


(3.24) 


This defines a new feedback system defined below: 


AX = BAX, 

where 

B = B + AB. 


(3.25) 

(3.26) 


The elements of B(t) are compared in Figures 3.1 and 3.2, where they are 
compared with the corresponding elements of B(t). 

The use of the polynomial gains represented by B(t) resulted in the 
end condition errors expressed by equation (2.17). Similarly, the follow- 
ing expression is obtained for B(t): 


AR = U(to) AX(t 0 ) 


L n 


+ / F (t) Af/m dt + / F x (t) 5X dt. 


/ 


(3.27) 
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The reduction in the effect of Af/m on Ar and A0 as a result of 
changing the gains from B(t) to B(t) is vividly illustrated in Figure 
3.3 where the elements of F 0 (t) are compared with the corresponding 
elements of F 0 (t). A comparison of the integrals of the absolute values 
reflects the following reduction in upper bounds. 



— ~ 


— ” 

1 ^ 1 
Ifol 

dt = 

.666 km sec 2 /m 

0Q> 

O 


*.258 deg sec 2 /m 

1_ J 


— — 



— — 


— — 

l?ol 


4.135 km sec 2 /m 


dt = 


0Q I 
O 


.550 deg sec 2 /m 

l— J 


_ — 


(3.28) 


(3.29) 


The reduction in the effect of Af/m could be expected because this 
was the criterion employed to determine the coefficient changes defined 
in equations (3.24). However, these coefficient changes are also 
expected to alter the effects of AX^t-) and AX. A comparison of the 
effect of these variables by using B(t) instead of B(t) is reflected in 
the following comparisons as well as in Figure 3.4 where the elements of 
Fi(t) are illustrated. 


.017670 

.019729 

.001184 

.002432 

-.001861 

-.000221 

-.001149 

-.000349 


(3.30) 


u(t c ) = 

.045516 

.033588 

.017544 

.009570 


-.007110 

-.000262 

-.003038 

-.000534 


(3.31) 


r” 

I£il 


.311 km/deg 



dt = 


/ 

0 

l£J 


_. 04 


(3.32) 
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(3.33) 



If if 


.498 km/ deg 


dt = 


liil 


.12 


The weighting functions of this system could again be linearized 
and further iterations performed; however, the results already obtained 
are sufficient to demonstrate the effectiveness of this operation. The 
additional improvement afforded by the inclusion of thrust acceleration 
feedback will be investigated. 

The expression for Al defined in equation (3.27) was derived for 
the following function. 


AX = B (t) AC + SX. 


The effect of a feedback of Af/m would give the following expression: 


AX = B(t) AX + e(t) Af/m + &X. 


(3.34) 


The expression corresponding to equation (3.27) would then be 


AR = U(t Q ) AXq + j t - 0 
t 


r n r 

/ [F 0 (t) + e(t) Fi(t) ] Af/m dt + / Fi(t) &X dt. 

t„ 


(3.35) 


Choosing e(t) as a quadratic function of time gives 


e(t) = e D + ejT + e 2 T 2 . 


(3.36) 


27 











Time [sec] 



FIG. 3-4. EFFECT OF COEFFICIENT CHANGES 
ON WEIGHTING FUNCTIONS OF A x 
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eQj^ej. and e 2 ^can be determined to minimize the weighted sum of squares 
of f 0 + e(t) f]_ and g Q + e(t) g l5 simultaneously. This was carried out 
for this example with the following results. 


e Q = -14.2412''! 
ei = 9.8851 ) . 

e 2 = - 1 . 4328> 


(3.37) 


The effect of Af/m is determined by the elements of the vector 


F 0 (t) + e(t) F 2 (t) . 

The integral of the absolute value of these functions is not significantly 
smaller than that for F 0 as shown below. 


— — 


— — 

|f Q + e (t) f x | 


.543 km sec 2 /m 


dt = 


|g Q + e(t) gij 


.190 deg sec 2 /m 



l?ol 


.666 km sec 2 /m 


dt = 


51 


.258 deg sec 2 /m 



(3.38) 


Although the improvement reflected by the above comparison is slight, 
a comparison of the functions themselves which appear in Figure 3.5 
shows that each of the elements of Fo(t) + e(t) Fi(t) changes sign 
which tends to cancel the effect of any bias of Af/m. Since a bias in 
Af/m is a likely occurrence, this will be considered sufficient justifica- 
tion for including the feedback of Af/m. 
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IV. CANCELLATION OF THE EFFECT OF INITIAL CONDITIONS 


Although the effect of initial conditions on AR is not consider- 
able, the linear portion is known at second stage ignition. It is 
possible to define a time function to cancel this effect. Suppose the 
thrust angle should be changed by an amount AX Q (t). The expression for 
AR would then include the following two terms in addition to others. 


u n 


4R = U(t Q ) AX(t 0 ) + / F x (t) AX^t) dt 4- 

U 

t^ 


(4.1) 


Let AX (t) be defined as follows: 


Z3C 0 (t) = -P(t) M _1 (t 0 ) U(t Q ) £X(t 0 ). 


(4.2) 


Substitution into equation (4.1) gives 


L n 


AR = U(t Q ) AX(t Q ) - J F x (t) P(t) dt M" 1 (t Q ) U(t Q ) AX(t 0 ). (4.3) 

‘o 

Equation (2.3) gives the following definition: 


-A 


M(t ) = F i(t) P(t) dt. 


Substituting this into equation (4.3) gives 


AR = U(t Q ) AX(t 0 ) - U(t Q ) AX(t 0 ) + ... , 


(4.4) 


NOT 


f\LN\ED- 


preceding 


33 


and the effect of initial conditions represented by U(t Q ) AX(t 0 ) is 
cancelled. The problem is reduced to determining a time function Ak(t) 
to add to the commanded thrust angle deviation so that the actual thrust 
angle will be changed by an amount AXo(t) after including the effect of 
the feedback. This can be expressed as follows: 


£X 0 (t) = B(t) AX(£X Q ) + Ak(t), 


(4.5) 


Thus, Ac(t) can be written 


£k(t) = AX 0 (t) - B(t) AX(AKq) 


(4.6) 


The linear expression for AX(£X 0 ) is 


AX^) = / U(t,t x ) F 1 (t 1 ) AK^ta.) dt x . 


(4.7) 


Substituting this expression into equation (4.6) and including the 
expression for AX Q (t) defined in equation (4.3) gives the following 
expression for £k(t): 


£c(t) = -P(t) + B(t) / U(t,ti) F 1 (t 1 ) P(t x ) dt x M _1 (t ) U(t Q ) £X(t 0 ) 


(4.8) 


This can be expressed as 


£c(t) » K(t) AX(t ) , 


(4.9) 


where 


K(t) = -P(t) B(t) / U (t , t x ) FjXtj.) P (t ]_) dt x M ” 1 (t Q ) U(t o ).(4.10) 
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K(t) is a 1 x 4 matrix of time-variable elements which will be 
fitted by some polynomial expressions. It will be convenient to use 
the same degree polynomials as will be required to fit the standard 
X with the feedback. Consequently, the time function required to 
allow the state variables to be measured directly, rather than in 
terms of differences, will be determined and fitted first. K(t) will 
then be fitted to the same degree polynomials. 

Except for the correction for initial conditions, the deviation in 
X which is predicted by the guidance function is 


XX = BAX + e(t) Af/m. 
This can be expressed as 


(4.11) 


X = BX + e(t) f/m + 



BX - e (t) 

b 



(4.12) 


Equivalently, 


X = BX + e(t) f/m + k(t). 


where k(t) is a polynomial in time. Since B(t) is composed of elements 
which are quadratic time functions and the variables x and y are expected 
to be at least quadratics, it appears that a function of at least fourth 
degree might be required. Since it fits satisfactorily, a sixth-degree 
polynomial is used in this example: 


k(t) = ko + k X T + k 2 x 2 + k 3 T 3 + k 4 r 4 + k 5 T 5 + k 6 x 6 . 


(4.13) 


The following numerical values were determined 

k Q = -2354.749911 > 
k x = -2267.341539 
k 2 = -738.009535 
k 3 = -39.95716990 > . 
k 4 = 5.34562706 
k 5 = -1.09387216 
k 6 = .10090398 - 


(4.14) 
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These coefficients were determined by the method of least squares 
with slight adjustments to ensure that AR = 0 for the standard trajec- 
tory. 


The effects of the initial conditions, similarly, were cancelled 
by fitting the four time functions of K(t) defined in equation (4.9) 
to sixth-degree polynomials of time. Slight adjustments in these 
coefficients were also made to ensure the following condition: 


J Fi(t) K(t) dt = -U(t Q ), 


(4.15) 


where K(t) is composed of functions fitted to approximate the elements 
of K(t). The time function £(t) can then be defined as follows: 


k(t) = k(t) + z£(t), 


(4.16) 


where 


£k(t) = Ak 0 + £kx T 4- £k 2 x 2 + 44 c 3 t 3 + £k 4 T 4 + Z4 c 5 t 5 + Z^k 6 T 6 . (4.17) 


These coefficients can be determined at second stage ignition by 
the following matrix operation. 


^o 
4^k x 

4k 3 

£k 4 

^6 


4K AX(t c ). 


(4.18) 


The values of the elements of ZJC which have been determined for 
this application are shown in the next section where a complete defini- 
tion of the initial and final conditions and the guidance function are 
listed with the results obtained from actual integration. 
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V. RESULTS AND CONCLUSIONS 


The guidance function derived in the preceding sections resulted 
in the following expressions: 

X = B(t) X(t) + e(t) f/m + k (t) , (5.1) 

B(t) = [a(t) S(t) £(t) d(t)] 

a(t) = a Q 4- a x T + a 2 x 2 

b(t) = b Q + b lT + b 2 X 2 

\ A - /S , A p 

c(t) = C 0 + C X T + C 2 T^ 

A A /N A r-j 

d(t) - d D + dj^T + d 2 x 2 

e(t) = e D + e^x + e 2 x 2 

k(t) = K(t) + Ak(t), 

k (t) = ko + kit + k 2 x 2 + k3t 3 + k 4 x 4 + ksX 5 + ks^ 6 
Ak(t) = £k 0 + £k x x + Ak 2 x 2 + Ak 3 x 3 + A 4 x 4 + AksX 5 + Ak 6 T 6 , 


where 


and 


x = 


t - t 
o 


100 


Ak 0 

Akj^ 

Ak 2 

Ak 3 

Z ^ 4 

Aks 

Ak6 


- AK AX(t 0 ), 


X(t 0 ) 


153.983 km 
6435.878 km 
2818.329 m/sec 
988.358 m/sec 
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The numerical values obtained for this example are as follows; 


/s 

a o = 

.104259 

a l = 

-.168328 

/N 

a 2 

= 

.122479 " 

/N 

^0 = 

.37021 

II 

H 

<rQ 

.315723 

b 2 

= 

.121797 

/s 

c o = 

.01551 

/N 

C 1 “ 

.041162 

c 2 

= 

-.003999 ) 

II 

o 

< ? TD 

.062320 

/N 

d x = 

.067045 

A 

= 

.000353 1 

e o ” 

-14.2412 

e l = 

9.8851 

e 2 

= 

-1.4328 J 


k Q = -2354.749911 
k x = -2267.341539 
k 2 = - 738.009535 
k 3 = - 39.957169900 ) 

k 4 = 5.34562706 

k 5 = -1.09387216 

k 6 = . 10090398 j 


(5.3) 


-.160198 

.249259 

-.246918 

-.043367 

-.072324 

.150760 

-.132564 

-.021482 

.229492 

.008411 

.152094 

.043929 

-.369459 

.224627 

-.376247 

-.082579 

.202657 

-.075912 

.180285 

.042941 

-.048929 

.021807 

-.045447 

-.010541 

.006095 

- . 002447 

.005513 

.001300 


(5.4) 


The above coefficients were derived for use with x and y measured 
in kilometers, x and y measured in meters per second, and f/m measured 
in meters per second 2 . In addition, they were derived for the follow- 
ing system of differential equations: 


x = f/m sin X + x 
g 


y 


f/m cos X + 




where 


x 

x g = 7 *• 


y g = J g. 


s o r S 


g g 0 = 9.81 m/sec 2 , r Q « 6,370 km. 


f/m = 


8.78065 


1.3751 - . 20888t ’ 


T = 


t - t 
C 

100 


The following initial conditions existed for the standard trajectory: 


x Q = 153.98343 km 

y Q =6435.8783 km 

x =2818.3294 m/sec 
o 

y Q = 988.35767 m/sec 

t = 146.815 sec. 
o 
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Cutoff was assumed to occur at the following velocity: 


V = 7792.5746 m/sec, 
c 

at which time the following values should exist for r c and 0 C . 


r = 6555.200 km 
c 

e c = 90 °. 


Table 5.1 presents a sample of twelve different second stage 
initial conditions caused by variations in first stage performance. 

The initial conditions are measured as deviations from the standard. 
Tables 5.2 and 5.3 show the resulting deviations in Ar and AG, respec- 
tively. The burning time of these examples were compared with those 
required for the corresponding calculus of variations solutions. 

Table 5.4 shows the additional propellant required for this guidance 
function beyond that which was required by the calculus of variations 
solution. Negative values are readily explained by the fact that 
errors in end conditions resulted in a different mission. The effect 
on burning time of errors in end conditions is readily determined from 
the following equation derived in Reference 2. 


tn 


At = -AjAr - AjAG + / [A2 + A3JE2 + A 2 g2] AX.^ 


(5.5) 


where 


Ai = -.08734 sec/km 
A 2 = .8589 sec/deg. 

For Example No. 13, with +1% F, W second stage perturbation, the follow 
ing errors were obtained: 


= -.161 km 
AG = .013 deg. 
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From this, the following variation in burning time should occur. (Only 
first order terms are considered.) 


At = -(-.08734) (-.161) - .8589(.013) = -.0252 sec. 


The flow rate for this case was W = 210 lb/ sec. Thus, in order to 
meet the end conditions indicated, additional propellant, ZW, should 
have been required. In this case, 


ZW = -5.3 lbs. 


The end conditions attained by the guidance function should have 
required 5.3 lbs of propellant less than that required by the calculus 
of variations solution to reach the nominal end conditions. However, 
Table 5.3 shows that the guidance function saved only 4 lbs of propel- 
lant, a net increase of 1.3 lbs more than the calculus of variations 
solution would have required to fulfill the same end conditions 
actually met by the guidance function. However, since the actual con- 
cern is the amount of fuel required to reach cutoff conditions from 
whatever cause, the comparison of fuel was made with respect to cal- 
culus of variations solution to nominal end conditions. 

In conclusion, it should be pointed out that the purpose of this 
report is not to propose a guidance function but rather to demonstrate 
several numerical methods on a rather intricate problem. These methods 
are applicable to a wide variety of problems, and it is felt that the 
results of this particular application demonstrate their effectiveness. 
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TABLE 5.1 

Deviations in Initial Conditions 


Example 

No. 

First Stage 
Deviations 

Ax (km) 
o 

Ay Q (km) 

Ak Q m/sec 

Ay o m/ sec 

1 

None 

0 

0 

0 

0 

2 

+5000 lbs 

- .93152 

-1.1623 

-30.2692 

-22.78749 

3 

-5000 lbs 

. 94444 

1.1824 

30.8335 

23.36393 

4 

Engine #2 
out at 100 sec 

17.28599 

5.5463 

- 4.6494 

-58.34471 

5 

Tail Wind 

.10183 

1.8180 

2.7670 

25.69483 

6 

Head Wind 

-.34235 

- .4055 

-5.3173 

- 5.40216 

7 

Left Cross Wind 

-.10434 

- .1788 

-1.5250 

- 2.33100 

8 

Right Cross Wind 

.09092 

.2155 

j 

1.4240 

2.87015 

9 

-1% W 

4.76116 

1.1588 

48.5244 

- 4.54684 

10 

+1% W 

-3.31692 

- .6959 

-39.7298 

3.74070 

11 

+1% F 

1.19032 

1.9245 

29.8569 

29.73483 

12 

-17. F 

-1.22157 

-1.9050 

-30.7597 

-29.22741 

13 

This example did not involve any first stage deviations. It 
is the case where the actual thrust angle exceeded the angle 
predicted by the guidance function by +1/2°. (SX = 4* 1/2°). 
Since it is combined with other second stage perturbations, 
it is listed with the first stage deviations. 






























TABLE 5.3 


A0 (degrees) 


PERTURBATIONS 
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